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Energy density of the Glasma 
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The initial energy density produced in an ultrarelativistic heavy ion colhsion can, in the color 
glass condensate framework, be factorized into a product of the integrated gluon distributions of the 
nuclei. Although this energy density is well defined without any infrared cutoff besides the saturation 
scale, it is apparently logarithmically ultraviolet divergent. We argue that this divergence is not 
physically meaningful and does not affect the behavior of the system at any finite proper time. 

PACS numbers; 24.85.+p,25.75.-q,12.38.Mh 



I. INTRODUCTION 

The matter produced at central rapidities in a heavy 
ion collision is dominated by the small x partons in 
the wave function of the high energy nuclei. These de- 
grees of freedom can, because of their high occupation 
numbers, be described as a classical Weizsacker- Williams 
color field. The source for this field is formed by the large 
X partons, which are seen by the small x ones as classi- 
cal color charges. The nonlinear interactions between 
the small x gluons give rise to gluon saturation, and the 
wavefunction is described by an energy (or x) dependent 
saturation scale. This way of understanding the small x 
wavefunction is known as the color glass condensate. A 
model incorporating these physical features was written 
down by McLerran and Venugopalan (MV) [11113. 

The initial transverse energy and gluon multiplicity in 
a collision of two sheets of color glass in the MV model 
has been calculated to all orders in the gluon field al- 
ready some time ago 0, IE IS S B I3- Recently there 
has been a renewed interest in the very early time be- 
havior of these classical "Glasma" gluon fields [Tol ll]| , 
in the context of pair production from the classical back- 
ground field |ll|il|ll|ll|ll|ll|ll|l|and parity vi- 
olation through the Chern-Simons charge density of the 
fields [20I |2]| . More attention has also been paid to the 
3-dimensional energy density (instead of the energy per 
unit rapidity) of these field configurations as a quantity 
that could be directly related to the initial conditions of 
hydrodynamical calculations [2^ |2^ . 

The purpose of this note is to clarify some properties 
of the initial energy density of the gauge fields in the MV 
model. We shall first, in Sec. ^ demonstrate that the 
initial energy density can be completely factorized into 
the product of the gluon distribution functions of the 
colliding nuclei. Going to finite proper times, or looking 
at the multiplicity, will change this factorization into a 
convolution of the unintegrated gluon distributions. We 
shall then, in Sec. IIIII go on to discuss the known prop- 
erties of the correlator of two pure gauge fields involved 
in the initial energy density and, in Sec. IIVI try to un- 
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FIG. 1; The color fields in different regions of spacetime. In 
regions (I) and (2), where only one of the nuclei has passed 
by, the field is the pure gauge field of this one nucleus. Inside 
region (3) the field is known numerically, but the initial con- 
dition for the field on the light cone r = is determined by 
the two pure gauge fields. 

derstand the behavior of large pt modes with the help 
of the lowest order perturbative solution of the classical 
field equations. 

II. INITIAL CONDITION 

Because of their high speed and Lorentz time dilation, 
the large x degrees of freedom are seen by the low x fields 
as slowly evolving in light cone time. They form classical, 
static (in light cone time) sources on the light cones: 

J" = <5^+P(i)(xT,x-) + J'^-p(2)(xT,a;+), (1) 

where the support of the sources around the light cone 
must be understood as being very close to a delta func- 
tion: P(i^2)(xt, 2;*) ^ S{x^). We shall work here in 
the Schwinger gauge Ar = {x^A^ +x~A'^ ) /r = 0, in 
which the current Q is not rotated by the soft classical 
field; in more general gauges Eq. ^ should be dressed 
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FIG. 2: Components of the color field, computed on a 512'^ 
lattice with g^fiRA ~ 67.7. 



in covariant gauge and Ref. 'll*! for another formulation.) 
This initial condition has a simple expression in terms of 
the pure gauge fields of the two colliding nuclei © : 
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Note that the metric in the {T,ri,'x.T) coordinate system 



is = diag(l,-T^,-l,- 
the Schwinger gauge 



-1) so that An 



-T^A''. In 



gauge field are related by A^ — ix^A"^ 



the ± components of the 
^A"^. Because of the 
explicit time dependence in the metric A'' corresponds, at 
T = 0, to the z-component of the chromoelectric field. At 
the initial time the only nonzero components of the field 
strength tensor are the longitudinal electric and magnetic 
fields and consequently the energy density is given by 



by Wilson lines to maintain its covariant conservation. 
The Weizsacker- Williams fields describing the softer de- 
grees of freedom can then be computed from the classical 
Yang-Mills equation 



£(t = 0) = lim 



1 dE 



(2) 



In the light cone gauge the field of one nucleus is a pure 
gauge outside the fight cone (see Fig. ^) 

^il,2)i^T) - ^t/(i^2)(xT)a.C/(\ 2)(xt), (3) 

where the SU(3) matrices U(^jyij{xx) are determined from 
the color sources as 

t^(i,2)(xT) = -Pexpj-ig J dx=^:^p(i,2)(xT,a;=^)| . 

(4) 

In the original MV model, which we shall be using here, 
the color charge densities are stochastic Gaussian ran- 
dom variables on the transverse plane: p^j^-)(xT, a;^) = 

S{x-)p"'{xt) with 

{p^xt)p'' {yr)) - g^p^S^'d^^T - Yt), (5) 

where the density of color charges g^p is, up to a numeri- 
cal constant and a logarithmic uncertainty, related to the 
saturation scale Qs- 

The initial conditions for the fields in the future light 
cone between the two colliding sheets were derived and 
the equations of motion solved to lowest order in the fields 
in [2J, |2^ |2g (see also Ref. ^7} for the same calculation 



0+ T d^XT drj 
1 



-Tr F F- 
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(10) 



Let us introduce a shorter notation for the correla- 
tion function of the pure gauge field of the nucleus when 
averaged with the distribution We shall define the 
correlation function G{pt) by 

A'r^^{pT)A^;^\ciT)) ^ 

(2^)2^™„^ab^2(p^ + c^t)^G{pt). (11) 

The index in parentheses (m) refers to the two colliding 
nuclei, which are, naturally, independent of each other, 
thus the (5"*" in the correlator. The correlator must also 
be diagonal in the color index (5"'') because there is no 
preferred direction in color space present in the prob- 
lem. We are assuming translational invariance on the 
transverse plane ((5^(pT + qT)), which is justified because 
we are only interested in momentum scales much larger 
than the nuclear geometry effects which break this in- 
variance (meaning that we are assuming \pt\ S> 1/Ra)- 
The transverse spatial index structure piPj is the only 
one consistent with rotational invariance on the trans- 
verse plane (again, at momentum scales |pt| ^ 1/Ra 
there is no preferred direction in the system to break 
this invariance). 

Using the notation of Eq. Hll|) the two terms in the 
energy density l{Tn|l become 
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(12) 
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and the final result factorizes into 
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Equation H14|l is our main result. The initial energy den- 
sity factorizes completely into a product of two terms, 
both of which only depend on the properties of one sin- 
gle nucleus. This happens only strictly at r = 0. 

Note that due to rotational invariance the initial en- 
ergy density in the magnetic field, Eq. (|12(l and the elec- 
tric field, Eq. (|13|l are equal. The discretized version of 
the computation on a transverse lattice breaks this rota- 
tional invariance (for an explicit expression see the lattice 
perturbation theory result in Appendix B of Ref. 3). 
This violation is largest for the momentum modes near 
the edges of the Brillouin zone. As we shall discuss in the 
following, for larger proper times these modes do not af- 
fect the energy density any more, and the energy densities 
in the longitudinal electric and magnetic fields approach 
each other, as can be seen in Fig. [2| 



III. PROPERTIES OF THE GAUGE FIELD 
CORRELATOR 

Let us then recall some known properties of G(pt), 
the correlation function of the pure gauge fields defined 
in Eq. Hll() . In light cone quantization it is related to the 
unintegrated gluon distribution function^ 



rG(k7 



(15) 



Note, however, that our G(pt) is equivalent to the unin- 
tegrated gluon distribution used to compute gluon pro- 
duction in pA-coUisions only in the weak field limit. We 
refer the reader to Refs. [23, 113, |^ for a discussion of 
the difference. Our G(pt) is equivalent, up to the nor- 
mahzation, to (j)^"^ of Ref. (23]. 

The correlator Gfpr) has been analytically evaluated 
in several papers [sj Ej, Q . The result is expressed in 
closed form in coordinate space as 



G(xt) 



5'A^cX^ 



1 
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(16) 



^ The notation and numerical constants at this point are very con- 
fusing. Here an attempt is made to follow Ref. 28] (in particular 
Sec. 2.4) and Ref. Q. 



Here A is an infrared cutoff that one must introduce in 
order to invert the 2 dimensional Laplace operator. The 
same function can also be measured in the numerical 
setup used to compute the glasma fields 0, IE IE S ■ 
The numerical procedure used is not exactly equivalent 
to the calculation leading to Eq. (|16|) , because the source 
in Eq. Q is taken as exactly a delta function on the light 
cone and the infrared divergence in inverting the Laplace 
operator is effectively regulated by the size of the 
lattice. These differences can, however, be absorbed into 
the infrared cutoff A, and the numerical evaluation (see 
in particular Fig. 3 of Ref. jlEl) of the correlator agrees 
with the behavior of Eq. H16|l . Note that to derive the 
correct initial conditions, Eqs. (|3J), © and ([TJl, it is 
essential to consider the source as spread out in the lon- 
gitudinal coordinate. Only when this is done can one, in 
practice, take the source as a delta function on the light 
cone when evaluating the Wilson line, Eq. 

Let us then estimate the behavior of the correlator 
G(pt) in momentum space. For small momenta G(pt) 
diverges, but the divergence is only logarithmic and thus 
integrable. This is the essential feature of gluon satu- 
ration; bulk quantities that are sensitive to the harder 
modes in the spectrum, such as the energy density, are 
infrared finite when the nonlinear interactions are taken 
into account fully. 

For large momenta G(pt) has a perturbative 
tail behaving as I/pt'^, meaning that the integral 
/ d^PTG{pT) and thus the initial energy density are 
seemingly ultraviolet divergent. This can be seen 
equivalently as the logarithmic ultraviolet divergence in 
Eq. 



(27r)2 



G(p, 



lim G(xt) 



27r5^ 



{g /i) lim In 



xrHO 



|xt|A- 



(17) 



We must emphasize that although Eq. (|17|l involves, for 
dimensional reasons, the infrared cutoff A, it corresponds 
to a divergence from large transverse momentum, or 
small distance, modes. 

The initial energy density of the glasma is infrared fi- 
nite, but seemingly logarithmically ultraviolet divergent. 
There are two reasons why this divergence is fundamen- 
tally not a problem for the physical picture of the glasma. 
The first reason is that, as can be seen from Eq. (|15|l . the 
divergence corresponds to large and thus, for a fixed 
energy, large x modes in the wavefunction. These are 
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FIG. 3: The energy density for different lattice spacings with 
^iRa = 67.7 The continuum limit corresponds to g^a 0. 
The fields configurations have been calculated on transverse 
lattices of f 00^ to 700^ points. The errors are statistical from 
between 20 configurations for the smaller lattices to 5 config- 
urations on the largest ones. 



degrees of freedom that were, by our initial assumptions, 
not meant to be included in the classical field in the first 
place. It would therefore be physically well motivated 
to regulate them with an ultraviolet cutoff Auv <^ Qs, 
and then match this cutoff with whatever way one treats 
these hard collisions. This is indeed the approach advo- 
cated e.g. in Ref. The other reason for not worrying 
about the ultraviolet divergence is that, as we shall ar- 
gue in the following, the energy density of the system at 
later proper times r 1/Qs has a finite Ayv — *■ oo limit. 
Thus if one regulates the ultraviolet divergence in any 
convenient way and proceeds to solve the equations for- 
ward in time, the cutoff no longer significantly infiuences 
the later time evolution of the glasma fields. 



IV. PERTURBATIVE COMPARISON AND THE 
APPARENT UV DIVERGENCE 



FIG. 4: The energy per unit rapidity /e = -^^^j-^^ 

(energy density divided by proper time) for different lattice 
spacings. The data are from the same numerical calculation 
as Fig. 
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FIG. 5: The energy density at proper times r of 0.2, 0.5 and 
1.0 in units of 1/g^fJ. as a function of the lattice spacing a, 
also in units of 1/g^fi. The data are from the same calcula- 
tion as Fig. 13 One sees that the energy density converges to 
the continuum limit {g^fJ.a 0) better at later times. The 
straight lines are linear extrapolations. 



In the perturbative (lowest order in the source charge 
densities) solution |24Ll25Ll26ll27l| the field amplitudes be- 
have, in the two dimensional Coulomb gauge, like Bessel 



m) is 

functions ^ 2 at f at 2 t\ f d^ky d^pT„. s 

e = -gN^ (TV, 1) y ^^G{pT)Gikr) 

A{T,kT) = MT = 0,kT)M\^T\r) (18) X [j2(|pr+kT|T) + j2(|pT + kT|T)] . (20) 

2t 

^7]{tt^t) = i^^i''' = Oi ki-) Ji (Ikj-lr). (19) Xo lowest order in the sources the pure gauge field cor- 
relator is 

The energy density corresponding to this perturbative , 1 (.g^/i)^ 

solution (this time dependence is also derived in Refs. [3^ 3^ Pt^ 
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and, using the asymptotic behavior of the Bessel func- 
tions in Eq. (|20|l . the energy density for large times re- 
duces to 

c(.^»)4/d»k,-j-^|M (22) 

with the Bertsch-Gunion [s^ type multiplicity resulting 
from the lowest order solution plIMlM l27| : 

dN _ ttR\ iVe(A^c^- l)gV^ 1 
dyd2kT ~ (27r)2 ^ 

f d^pT 1 

J (2^)2 p^2(k^_p^)2- ^ ^ 

Also in the full nonperturbative solution the most ul- 
traviolet modes behave in the same way 0, IE US ■ Let 
us now assume that we have computed the initial en- 
ergy density with some ultraviolet cutoff Auv (such as 
the inverse lattice spacing in the numerical calculation). 
The initial energy density Eq. (|14l) depends on this cut- 
off as In^Auv, as can be seen in the numerical result 
in Fig. 13 After a time r > 1/Auvi however, the time 
dependence of the modes near the cutoff changes from 
the initial Ai(r) = Aj(0) (l + 0{t'^)) to the asymptotic 
regime Ai ^ l/|kT|T and the contribution of these modes 
to the total energy density is suppressed by an additional 
power of the momentum, making the energy density fi- 
nite in the limit Ayy — >■ oo. After a time r > 1/Auvi 
or T > a in the lattice computation, the energy density 
then converges to a value that is independent of the lat- 
tice spacing, as shown in Fig. ^ The larger the proper 
time that one is looking at, the better the convergence 
in the continuum limit is. Figure [S] shows the continuum 
extrapolations of the energy densities at g'^fiT = 0.2, 0.5 
and 1.0. 

By turning this argument around one can see that if 
one takes first the continuum limit (Auv — > oo or a ^ 0) 
for a finite r, the continuum extrapolated energy den- 
sity will behave as In^ 1/r. Thus if the limits are taken 
in this order, the energy density is indeed finite for all 
T > 0, but diverges logarithmically at r = 0. Note that 
this divergence is so weak that the energy per unit ra- 
pidity (re) is still zero at r = 0. Thus we se that the 
solution of the field equations is well defined as an ini- 
tial value problem only in the presence of an ultraviolet 
cutoff in the transverse momenta, but for later times the 
system loses memory of this cutoff. Incidentally, because 
of this feature one can argue that introducing of a finite 
initial time (such as done in Refs. ^ to avoid the 
singularity resulting from broken boost invariance of the 
field configurations) does not really add another physical 
parameter into the model. 

V. ENERGY DENSITY IN PHYSICAL UNITS 
AND CONCLUSION 

For concreteness let us finally try to express the re- 
sults shown in Figs. |31 El and |S1 in physical units. Due 



to the difficulty in fixing exactly the right value of the 
color charge density parameter g^fj, this is not necessarily 
straightforward. For RHIC energies one can argue, based 
on both the gluon multiplicity found in the numerical 
calculation i and counting the number of large x 
degrees of sources in the wavefunction from conventional 
parton distribution functions 26] that the relevant value 
would be g^fj. = 2 GeV. Other estimates, e.g. jl^], give a 
smaller value, so g"^!! = 2 GeV should be considered an 
upper bound for RHIC. For LHC energies the estimate 
based on parton distributions 0| gives g^fi ~ 4 GeV, but 
this is most certainly an overestimate, since the calcula- 
tion in Ref. |26J is based on parton distributions in the 
proton and does not take into account shadowing correc- 
tions. Another way of estimating the color charge den- 
sity is based on the small x scaling properties observed 
in deep inelastic scattering data [ia, El, Q EE El • The 
saturation scale in this scaling can then be related to 
the MV model color charge density 47, 48]. This line 

of thought leads to a scaling ^ {g^fJ-)^ ^ V^^, where 
A = A/(l + A/2) and_a fit to the HERA data 42] gives 
A ^ 0.288 and thus A w 0.25. The resuh for the color 
charge density at the LHC would be g^ii ~ 3 GeV, which 
is the value we will use in the following. 

As we have seen, the energy density strictly at t = 
is not the best quantity to look at. Let us instead 
estimate the energy density at the time r = 1/g^^. This 
is when, as can be seen from Fig. 0] the 1/r-decrease 
of the energy density seems to start. The simple linear 
continuum extrapolation of Fig. [Slyields e(T = 1/g'^fj,) = 
0.26(52^)4/(^2^ Using g = 2 we then get the estimates 
e{T = 0.1 fm) w 130 GeV/ fm^ for RHIC and e(r = 
0.07 fm) « 700 GeV/fm^ for the LHC. This estimate 
agrees with the values given in Ref. ^E] for t — 3/g^n 
when the ~ 1/r behavior of the energy density is taken 
into account. The uncertainty due to the unprecise value 
of g2/i in these numbers is quite large because of the 
power law dependence e ^ (g^^)'^. 

In conclusion, we have shown that the initial 3 dimen- 
sional energy density of the Glasma fields in the MV 
model can be expressed as a product of of the (integrated) 
gluon distribution functions of the colliding nuclei. Only 
the energy density at later times and the multiplicity 
involve convolutions of the unintegrated gluon distribu- 
tions, probing the kT-distributions in the wavefunctions 
of the colliding nuclei in more detail. We have recalled 
the known properties if the pure gauge field correlator 
Eq. Hll|) appearing in the initial energy density. As ex- 
pected from general gluon saturation arguments, the ini- 
tial energy density is infrared finite when the gluon fields 
are solved to all orders in the source. The energy den- 
sity strictly at t = is, however, ultraviolet divergent 
in the MV model. We show, both by a direct numerical 
calculation and by examining the time dependence of the 
ultraviolet modes, that this divergence does not persist 
when the equations of motion are solved to times larger 
than the inverse ultraviolet cutoff. 



6 



Acknowledgments 

The author would Uke to thank K. Kajantie, L. McLer- 
ran and R. Pries for discussions that led to writing this 



paper and R. Vcnugopalan for urging to actually write it 
up and comments on the manuscript. This manuscript 
has been authorized under Contract No. DE-AC02- 
98CH10886 with the U.S. Department of Energy. 



[1] L. D. McLorran and R. Venugopalan, Phys. Rev. D49, 

2233 (1994), [arXiv:liep-ph/9309289]. 
[2] L. D. McLerran and R. Venugopalan, Phys. Rev. D49, 

3352 (1994), [arXiv:hep-ph/9311205]. 
[3] L. D. McLerran and R. Venugopalan, Phys. Rev. D50, 

2225 (1994), [arXiv:hep-ph/9402335]. 
[4] A. Krasnitz and R. Venugopalan, Nucl. Phys. B557, 237 

( 1999) , [arXiv:hop-ph/9809433] . 
[5] A. Krasnitz and R. Vcnugopalan, Phys. Rev. Lett. 84, 

4309 (2000), [arXiv:hcp-ph/9909203]. 
[6] A. Krasnitz and R. Venugopalan, Phys. Rev. Lett. 86, 

1717 (2001), [arXiv:hep-ph/0007108]. 
[7] A. Krasnitz, Y. Nara and R. Venugopalan, Phys. Rev. 

Lett. 87, 192302 (2001), [arXiv:hep-ph/0108092]. 
[8] T. Lappi, Phys. Rev. C67, 054903 (2003), [arXivihep- 

ph/0303076]. 

[9] T. Lappi, Phys. Rev. C70, 054905 (2004), [arXivihep- 
ph/0409328]. 

[10] T. Lappi and L. McLerran, Nucl. Phys. A772, 200 

(2006), [arXiv:hcp-ph/0602189]. 
[11] R. J. Fries, J. I. Kapusta and Y. Li, arXivinucl- 
th/0604054. 

[12] F. GeUs and R. Venugopalan, Phys. Rev. D69, 014019 

(2004) , [arXiv:hep-ph/0310090]. 

[13] F. Gclis, K. Kajantie and T. Lappi, Phys. Rev. C71, 

024904 (2005), [arXiv:hcp-ph/0409058]. 
[14] T. Lappi, arXivdicp-ph/0606090. 

[15] F. Gelis, K. Kajantie and T. Lappi, Phys. Rev. Lett. 96, 

032304 (2006), [arXiv:hep-ph/0508229]. 
[16] H. Fuju, F. Gelis and R. Venugopalan, Eur. Phys. J. 

C43, 139 (2005), [arXiv:hep-ph/0502204]. 
[17] H. Fujii, F. Gclis and R. Vcnugopalan, Phys. Rev. Lett. 

95, 162002 (2005), [arXivdicp-ph/0504047]. 
[18] D. Kharzeev and K. Tucliin, Nucl. Phys. A753, 316 

(2005) , [arXiv:hep-ph/0501234]. 

[19] D. Kharzeev, E. Levin and K. Tuchin, arXivihep- 
ph/0602063. 

[20] D. Kharzeev, A. Krasnitz and R. Venugopalan, Phys. 

Lett. B545, 298 (2002), [arXiv:hcp-ph/0109253]. 
[21] D. Kharzeev, Phys. Lett. B633, 260 (2006), [arXiv:hep- 

ph/0406125]. 

[22] T. Hirano and Y. Nara, Nucl. Phys. A743, 305 (2004), 

[arXiv:nucl-th /0404039] . 
[23] T. Hirano and Y. Nara, J. Phys. G30, S1139 (2004), 

[arXiv:nucl-th/0403029] . 
[24] A. Kovncr, L. D. McLerran and H. Weigert, Phys. Rev. 

D52, 6231 (1995), [arXiv:hep-ph/9502289]. 
[25] A. Kovner, L. D. McLerran and H. Weigert, Phys. Rev. 

D52, 3809 (1995), [arXiv:hep-ph/9505320]. 



[26] M. Gyulassy and L. D. McLerran, Phys. Rev. C56, 2219 

(1997), [arXiv:nucl-th/9704034]. 
[27] Y. V. Kovchegov and D. H. Rischke, Phys. Rev. C56, 

1084 (1997), [arXiv:hep-ph/9704201]. 
[28] E. lancu, A. Leonidov and L. D. McLerran, Nucl. Phys. 

A692, 583 (2001), [arXiv:hep-ph/0011241]. 
[29] D. Kharzeev, Y. V. Kovchegov and K. Tuchin, Phys. 

Rev. D68, 094013 (2003), [arXivdicp-ph/0307037]. 
[30] J. P. Blaizot, F. Gelis and R. Vcnugopalan, Nucl. Phys. 

A743, 13 (2004), [arXiv:hep-ph/0402256]. 
[31] F. Gelis, A. M. Stasto and R. Venugopalan, arXivihep- 

ph/0605087. 

[32] J. Jalilian-Marian, A. Kovncr, L. D. McLerran and 
H. Weigert, Phys. Rev. D55, 5414 (1997), [arXivdicp- 
ph/9606337]. 

[33] Y. V. Kovchegov, Phys. Rev. D54, 5463 (1996), 

[arXiv:hep-ph/9605446] . 
[34] Y. V. Kovchegov, Phys. Rev. D55, 5445 (1997), 

[arXiv:hep-ph/9701229]. 
[35] A. Krasnitz, Y. Nara and R. Vcnugopalan, Nucl. Phys. 

A717, 268 (2003), [arXiv:hep-ph/0209269]. 
[36] Y. V. Kovchegov, Nucl. Phys. A762, 298 (2005), 

[arXiv:hep-ph/0503038] . 
[37] Y. V. Kovchegov, Nucl. Phys. A764, 476 (2006), 

[arXiv:hcp-ph/0507134]. 
[38] J. F. Gunion and G. Bcrtsch, Phys. Rev. D25, 746 

(1982). 

[39] A. Krasnitz, Y. Nara and R. Venugopalan, Nucl. Phys. 

A727, 427 (2003), [arXiv:hep-ph/0305112]. 
[40] P. Romatschke and R. Venugopalan, Phys. Rev. Lett. 

96, 062302 (2006), [arXiv:hep-ph/0510121]. 
[41] P. Romatschke and R. Vcnugopalan, Phys. Rev. D74, 

045011 (2006), [arXiv:hcp-ph/0605045]. 
[42] K. Golcc-Bicrnat and M. Wusthoff, Phys. Rev. D59, 

014017 (1999), [arXiv:hep-ph/9807513]. 
[43] K. Golec-Biernat and M. Wusthoff, Phys. Rev. D60, 

114023 (1999), [arXiv:hep-ph/9903358]. 
[44] A. M. Stasto, K. Golcc-Bicrnat and J. Kwiccinski, Phys. 

Rev. Lett. 86, 596 (2001), [arXiv:licp-ph/0007192]. 
[45] A. Freund, K. Rummukainen, H. Weigert and A. Schafer, 

Phys. Rev. Lett. 90, 222002 (2003), [arXiv:hep- 

ph/0210139]. 

[46] D. Kharzeev, E. Levin and M. Nardi, Nucl. Phys. A747, 

609 (2005), [arXiv:hcp-ph/0408050]. 
[47] E. lancu and R. Vcnugopalan, arXiv:hep-ph/0303204. 
[48] H. Weigert, Prog. Part. Nucl. Phys. 55, 461 (2005), 

[arXiv:hep-ph/0501087]. 



